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ABSTRACT 

Whether among the myriad tiny proto-galaxies there exists a population with similarities to present 
day galaxies is an open question. We show, using BlueTides, the first hydrodynamic simulation large 
enough to resolve the relevant scales, that the first massive galaxies to form are predicted to have 
extensive rotationally-supported disks. Although their morphology resembles in some ways Milky- 
way types seen at much lower redshifts, these high-redshift galaxies are smaller, denser, and richer 
in gas than their low redshift counterparts. From a kinematic analysis of a statistical sample of 216 
galaxies at redshift z = 8 — 10 we have found that disk galaxies make up 70% of the population of 
galaxies with stellar mass IO^^Mq or greater. Cold Dark Matter cosmology therefore makes specific 
predictions for the population of large galaxies 500 million years after the Big Bang. We argue that 
wide-field satellite telescopes (e.g. WFIRST) will in the near future discover these first massive disk 
galaxies. The simplicity of their structure and formation history should make possible new tests of 
cosmology. 


1. INTRODUCTION 

The search for “primeval” galaxies and the understand¬ 
ing that will come from studying them has long been a 
guiding principle for much of extra-galactic astrophysics 
and cosmological research. Theories for the formation of 
the first galaxies s pan early models o f monolithic collapse 
of massive bodies(Eggen et al.||1962|) to merger scenarios 
for the hierarchical formation of all structures!Press & 


Schechter 1974). Disk galaxies are arguably those that 
have generated the mo st interest, and were the focus of 
early theoretical work(Fall & Efstathiou 1980|). In the 
simplest galaxies, the structure of rotating disks should 
be directly linked to the angular momentum of the sur¬ 
rounding dark matter halo. Einding galaxies for which 
this is true, in the earliest environments possible would 
enable powerful tests of the dynamics, gaseous and ra¬ 
diative processes involved in galaxy formation. 

In current cosmological models, galaxies form from the 
gravitational co llapse of small perturb ations in the mat¬ 
ter distribution! Silk & Mamon 2012). This process in¬ 
volves both a hierarchy of merging structures and smooth 
accretion, so that early galaxies are predicte d to be 
morphologically irregular, clump y, and compact!Shimizu 
et al. 2014 Cen & Kimm 2014). This is supported by 
recent observational d ata on samples of galaxies at red- 
shift 2 ; = 8 and beyond!|Oesch et al.|20I0||Ono et al. 12013' 


eyondljj 

.f20I4| ) 


Curtis-Lake et al. 120 141). The volumes accessible to these 


studies, both coniputational and observational are how¬ 
ever thousands of times smalle r than those that w ill be 
probed by upcoming telescopes! Spergel et al.||20I5 ). 
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The earliest that large rotating disks h ave been seen 


observationally is at re dshifts 2 : = 2 — 3 !Genzel et al 
2006 Glazebrook|20I3 ), and this has required deep spec¬ 


troscopy on I0m%lass telescopes. At higher redshifts, 
galaxies in published samples are v ery compact, with 


half light radii ^ 0.5 kpc at 2 : = 8!Oesch et al. 2010 


Holwerda et al. 2014). The fraction of morphologically 


disturbed galaxies is also known to be high in deep fields 
!Menanteau et al. 2006). The volumes probed at red- 
shifts 2 ; > 3 and particularly the redshifts at the current 
frontier, 2 ; = 8 — 10 are however very small, with a sky 
areas of II sq. arcmin accessible t o the Hubble Extreme 
Deep Eield!Illingworth et al.|20I3) for example. This has 
meant that probing the highest luminosity and largest 
galaxies, has not been possible so far. It is not surprising 
that observational studies of the structure and morpholo¬ 
gies of the highest 2 ; galaxies have not revealed large disk 
galaxies, whether they exist or not. 

Theoretically, the situation is similar - it is not known 
what the prevailing Gold Dark Matter model predicts 
for the most massive galaxies at redshifts 2 ; > 8. Nu¬ 
merical simulations used to make predictions have so 


far! 

Jaacks et al.| 

2012 

Gen & Kimm||20I4) been limited 

to volumes 4U-R 

JO times smaller than our present work. 


or else zoomed simulation s concentrating on individual 
galax ies ! |Pawlik et al. 12011 ) or overdensities! Yajima et al 


2014). intriguingly, a dominant fraction of the early 


forming (dwarf-mass) galaxie s simulated so far show ev¬ 


idence of disk-like structure! Romano-Dfaz et al. 2011, 
Pawlik et al. |2011|). It is of prime importance to reach 


the regime of large galaxies as the early formation of lu¬ 
minous or massive objects (such as first bright quasars, 
and most massive galaxy clusters) is one of the most po¬ 
tent tests of models of structure formation in hierarchical 
cosmology. In order to do so, and simulate a statistically 
useful sample of galaxies, it is necessary to model a vol¬ 
ume of order several hundred megaparsecs on a side. This 
must be done with high enough mass and force resolution 
to allow studies of structure and morphology (See e.g.. 
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Figure 1. A sample of disk galaxies selected from the BlueTides simulation at redshift z = 8. We show galaxies identified kinematically 
to be rotating disks. For each galaxy we show a face-on and side-on view. The colors of top two rows represent the stellar surface density 
coded by the stellar age (blue to red: 0 ~ SOOMyear), face-on and side-on. The colors of the bottom two rows represent the star-formation 
surface density (normalized to 0 ~ 1). The face-on images have included the effect of dust-extinction. The darker central regions of the 
disks are mostly caused by the consumption of gas to form bulge stars. There is also an obscuration effect from the higher metallicity (and 
therefore dust content) of the central regions. 


Covernato et al. 
regime, with the 


f 007|)). We have been able to reach this 
lueT'ides run, a cosmological hydrody¬ 


namic simulation which we have just completed, and we 
describe the first results below. 

2. THE BLUETIDES SIMULATION 

BlueTides was carried out using the S moothed Particl e 
Hydrodynamics c ode M P-Gadget (see Springel (2005); 
Di Matteo et al. (2012|) for earlier versions of the code) 
with 2 X 7040^ (0.69 tnllion) particles on 648,000 Cray 
XE compute cores of the Blue Waters system at the Na¬ 
tional Center for Supercomputing Applications (NCSA). 
The simulation evolved a cube of side length 400 comov¬ 
ing Mpc/h to redshift z = 8, and is the largest cosmo¬ 
logical hydrodynamic simulation yet carried o ut in terms 
of memory usage, by an order of magnitude! Di Matteo 


et al 
mode 


2012[ ). The supernova feedb ack in BlueTides is 
ed after Okamoto et al. (2010), with the wind feed- 


back parameter = 3.7. 'This model has been shown to 
reproduce several key pro perties of low redshift galaxies 
( Vogelsberger et ar]|2014 ). We also include the effect of 
super-massive black holes, primordial and metal cooling, 
H2 molecular cloud and a non-uniform ionization back¬ 


this is shown in Feng et al. (2015). In the present case, 
we are most interested in the most massive galaxies with 
at least 10^ stellar particles, which contain enough parti¬ 
cles to study their rotational support. The corresponding 
stellar mass is 0.7 x IO^^Mq or more. 

We have made images of all galaxies present in the sim¬ 
ulation, by adaptively smoothing the star particle distri¬ 
bution to reveal as much structure as possible given the 
resolution available. We have also used the star forma¬ 
tion rate (SFR) of the gas particles in the simulation to 
make other images. The SFR can be related to the rest- 
frame UV (150 nm) luminosity. The UV luminosity is 
sensitive to the presence of dust, and there is evidence 
that the highes t luminosity early galaxies are signifi - 
cantly obscured! Wilkins et al. pM3||Cen & Kimm|2014D. 
We ac count for this dust attenuation using a screening 
model! Joung et al. 2009). This simple model assumes 
that the dust attenuation Auv ^ galaxy pixel is pro¬ 
portional to the metal-mass density in that pixel. There 
is only one free parameter, which we normaliz e by set¬ 
ting Auv = 1 for observed Muv = — 21 galaxies! Wilkins 


ground: we refer to Feng et al.|!2015) for the parameters 
and other aspects of the simulation in more detail. We 
also show there that the global star formation rate pre¬ 
dicted is consistent with current observations at redshifts 
2 : > 8. The volume simulated is approximately 300 times 
larger than the lar gest observational survey at these red- 
shifts !fT4enti||2M^ . 

Galaxies were selected from the raw particle data (47 
TB per output time) using a friends-of-friends algorithm 
at a range of redshifts from 2 : = 8 — 13. The mass reso¬ 
lution of the simulation is 1.72 x IO^Mq per dark mat¬ 
ter particle, so that at redshift 2 : = 8, 36 million galax¬ 
ies were found containing 64 dark matter particles or 
more, above a mass threshold of IO^Mq. The BlueTides 
galaxy luminosity functions at the redshifts of interest 
here are fully consistent with current observations, and 


et al. 2013|). We have taken this approach, rather than 
performing a more detailed radiative transfer simulation 
because the current major un certainties lie in the dust 


content and attenuation curve! Wilkins et al.||2Q13 ). 


3. GALAXY MORPHOLOGY 

We turn directly to the morphology of the most mas¬ 
sive (Msteiiar ^ 0.7 X IO^^Mq) galaxics. We compute 
the inertia tensor of the star particles and use these to 
plot views of each galaxy parallel and perpendicular to 
the major axis. Studying them reveals that a significant 
fraction are visually disk-like, and surprisingly regular in 
shape. 

In Figure 1 we show images of 10 disk galaxies in the 
simulation at ^ = 8, both side-on and face-on for each. 
These are among those kinematically selected to be disks 
(see below). We can see that the disk components extend 
to a diameter of ^ 4 kpc, with a median half-light radii 
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Figure 2. Fraction of galaxies in the BlueTides simulation kine¬ 
matically classified as disks. Results are shown as a function of 
galaxy stellar mass, and at 3 different redshifts. 

^0.6 kpc for those shown. Although the 2 ; = 8 BlueTides 
disk galaxies have masses comparable to the Milky Way, 
their sizes are significantly smaller, by a factor of ^ 5. 
The disk galaxies also have mean gas fractions of ^ 70% 
at ^ = 8. The typical stellar mass/halo mass ratio of 
these systems is a few percent 5%, estimated from 
Friend-of-Friend halos). This is similar to systems at 
z=0, where the stellar to halo ratio peaks at 4% for halos 
in the same mass range (|Guo et al.||2010|). The star for¬ 
mation rates of these galaxies are ^ 100 MQ/year, orders 
of magnitude larger than th e main sequence galax ies of 
similar stellar mass at low-z (Whitaker et al. |201^ . The 
stellar populations in the central regions of each galaxy 
are significantly older than the outer disk, with a pro¬ 
portion of stars that formed as early as 2 : = 13, i.e. age 
300 Myr at 2 : = 8. 

These disk galaxies host supermassive blackholes, with 
blackhole masses ranging from a few times 10^ to a few 
times 10^ ^“^Mq. Interestingly, the host galaxy of the 
most massive blackhole (3 x 10 ^/i“^Mq) has a much more 
dominate spheroidal component {DjT = 0.29) than typ¬ 
ical galaxies with similar masses {DjT > 0.5). We plan 
to investigate the role of central blackholes in the forma¬ 
tion of disks in the near future. 

Kinematic decomposition of the galaxies lets us know 
if the disk-like structures seen in Figure 1 are rapidly ro¬ 
tating. We have used the standard technique! Governato 
et al. 2007) to determine the fraction of stars in each 
galaxy which are on planar circular orbits, and which 
are associated with a bulge. We note that this technique 
is not sensitive to the dumpiness of disks. At redshift 
2 ; = 8 we find that 70% of galaxies above a mass of 
10 ^^Mq are kinematically classified as disks (using the 
stand ard threshold of disk sta rs to total stars (D/T) ratio 
of 0.2(Governato et al. 2007). This percentage of disks 
changes with redshift as shown in Figure 2. In all cases it 
is much higher than the fraction of observed galaxies of 


Figure 3. Disk galaxy kinematics in BlueTides characterized by 
the ratio of circular velocity V to vertical velocity dispersion a. 
We show histograms of results for this V/a ratio computed using 
star particles in the disks of galaxies selected to be disks from their 
D/T ratio. Data on redshifts z = 8 and z = 9 are shown. 

comparable mass which are seen to be disk-like at ^ = 0, 
0.14 ± 0.03(|Buitrago et al.||2013). 

We note that the majority of gas-rich disks seen at 
high redshifts accessible so far t o spectroscopy (z=2-3) 
are cl u mpy, turbulent and t hick (Forster Schreiber et al. 
2009) (Genzel et al. (2006) is an exception), and it is 
natural to expect this to be the case with the BlueTides 
massive disks, which have high gas fractions. A relevant 
statistic we have computed is the ratio of circular veloc¬ 
ity, V to vertical velocity dispersion a of our disk galax¬ 
ies. We have computed this Vja ratio using information 
from the disk star particles and plot a histogram of the 
results in Figure 3. We can see that the values range from 
F/cr = 4 to 6, with a mean value of 4.98 ±0.04 measured 
from all redshifts. The galaxies are therefore definitely 
disks, having Vja well above the 0 — 1 expected for non¬ 
disk galaxies. This supports our conclusions made using 
the D/T ratio above. The Yja values are lower than 
the Vja ^ 10 expected for thin spiral galaxies seen at 
2 : = 0, and instead con sistent with gas rich turbulen t 
disks from earlier times (jForster Schreiber et al. 2009). 
Although this is reasonable, confirmation of this aspect 
awaits future higher resolution resimulations which are 
able to resolve the fine vertical structure of the disks. 

4. HALF-LIGHT RADIUS 


A number of recent studies! Oesch et al. 2010| |Ono 


__ _ il |Ono 

et al.|[2QT^ IHolwerda et ar]|2014[ [Kawamata et al.[r2014[ 

Shibuya et al.||2015p have begun to make measurements 

of galaxy structure and in particular sizes at the high¬ 
est redshifts. As the observational volumes are small, 
the magnitude range of observed samples varies from 
Muv = “21 to Muv = —18. We can therefore com¬ 
pare BlueTides galaxies with these measurements, bear¬ 
ing in mind that the brightest galaxies in the simula¬ 
tion are much rarer than can be observed in current sur- 
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Figure 4. Galaxy half-light radii at redshift z = S, 9, 10 in the 
BlueTides simulation (including disks and non-disks) and in the 
Hubble Space Telescope observations of Holwerda et al (2014), 
Kawamata et al (2014), Curtis-Lake (2014), and most recently 
Shibuya et al (2015). Colored Symbols: Galaxies in BlueTides Sim¬ 
ulation; We include galaxies with a stellar mass > 5 x 10® 

Upward wedges: Holwerda et al (2014) z > 9 samples; Downward 
wedges: Kawamata et al (2015) z = 8 samples. Circles: Curtis- 
Lake et al (2014) z = 8 undisturbed samples Muv > —20.5, with 
uncertainty estimated from their Figure 1. The half-light radius in 
BlueTides is measured from the face-on images. The right-axis is 
in units of the force resolution Cgravi = 0.24 Kpc at 2 ; = 8. 

veys. Using the standard observational algorithm SEX- 
tractoiQ we measure the half-light radius (ri/ 2 ) of each 
galaxy. This is shown in Figure 4, where we compare 
our redshift z = S results with observational data from 


([Kawamata et al. |2Q14), ([Holwerda et al.|2Q14 ) 
Lake et al. [2U14[ ) and (!5hibuya et al. ' |2U15[ ). 


(Curtis- 
We can 

see that the current observations cover the lower end 
of the magnitude range, and are in the same realm as 
(but not a good match to) our simulation results, with 
a mean ri /2 ^ 0.6 for Muv ^ —20.5 —21.5. Galaxies 

are therefore extremely compact at these high redshifts. 
The galaxy sizes in the simulation do appear to be on 
the low end compared to o bservations, somet hing which 
has also been seen at z=3(Joung et al. 2009) and could 
be a sig n of interesting physics, it has also been shown 


however (Curtis-Lake et al. 2014) that selection effects 


flux cuts and choice of size measurement algorithms can 
affect measured galaxy sizes by factors of 2 or more. 

The half light radius of the faintest samples 0.5Kpc) 
is only a factor of 2 of the gravitational smoothing length 
(0.24 Kpc at z = 8), which could also affect our results. 
In future work we will explore these effects in more de¬ 
tail with mock observations. Interestingly, the simulation 
predicts a positive size-luminosity correlation. 

Analytic models can be used to predict how disk sizes 
scale with redshift. For example, the simplest assump¬ 
tion is that the half-light radii scales with virial radius. 


which has the following redshift dependence 
Vyir OC OC (1 -h 


for constant halo circular velocity (e.g.,[Mo et al.[([1998[)). 
Using this relationship to extrapolate the dis k ^zes from 
the GEMS o bservational sample analyzed in SomerviH^ 


with 


)8| (in which the median for galaxies 
Mq at 2: = 1.05 is 2.1Kpc), this simplest model 


et al. (^2008| 

= lU^ 

predicts that disks of this mass should have Vd = 0.23 kpc 
at 2 ; = 8 and Vd = 0.17 kpc at 2 ; = 10. This is somewhat 
lower than disks of similar stellar mass in our simula¬ 
tion which have half-light radii about twice this value 
(see Figure 4 - these correspond to th e least luminous 
disks plotted). Somerville et al. (2008) have made more 
sophisticated rnodels which include evolution in the in¬ 
ternal structure of dark halos and adiabatic contraction, 
both of which mean that the trend of size evolution with 
redshift is expected to be weaker than in the simplest 


model s (also supported by observations, e.g. Shibuya 
et al. (2015)), and which therefore appear to describe 
our sminiated results better. 

5. CONCLUSIONS 

At redshifts ^ = 8 — 10 the Universe is about 1000 
times denser than the present day. In these very different 
physical conditions, and in the very short time available, 
massive disk galaxies are able to form in the cold dark 
matter model. They share many visual and kinematic 
characteristics with those that are seen when the Uni¬ 
verse was 20 times older. This means that an image one 
might have of the early Universe containing only small 
protogalaxies, merging clumps and irregular structures 
is not appropriate. 

An immediate question we can ask is how these disks 
did in fact form in the simulation. We have tracked parti¬ 
cles in the most massive galaxies (including those plotted 
in Figure 1) back in time to visually assess their forma¬ 
tion history. We have found that of the most massive 
20 galaxies at redshift 2 ; = 8, only one was the prod¬ 
uct of an obvious major merger between similarly sized 
progenitor galaxies. A popular pi cture of gas rich ma¬ 
jor m ergers producing some disks (Springel & Hernquist 
2005) therefore does not appear to be rele vant here. 

Eve n though t he Lagrangian volu mes ([Ohorbe et al. 
2014) and halos (Rossi et al. 2011[) at the rare peaks 


^ We use a Python rewrite of SEXtractor by [Barbary fc contrib-[ 
|utors| ( [2014| ). 


associated with these galaxies are likely more spherical 
the distribution of matter around galaxies is known to 
be highly anisotropic in the CDM model, the infall of 
gas cannot proceed quasi-s pherically as in the clas sical 
picture of galaxy formation(Fall & Efstathiou[[1980). In¬ 
stead, smooth infall and accretion build up the disks, 
with cold gas arriving along the dir ections of nearby fil- 
aments in the density distribution (Romano-Dfaz et al. 
2 ^. 

These massive disks will be fascinating objects to study 
wit h next generation te lescopes. The WFIRST satel- 
lite(Spergel et al.l 2015[) will have a field of view 200 
times that of the WL'C3 instrument on HST, making a 
planned deep (J=26.7) sky survey covering 2000 square 
degrees possible. Using the luminosity function in Blue¬ 
Tides and its evolution we predict that this survey should 
find ^ 8000 disk galaxies of the type we have kinemati¬ 
cally identified as massive disks, about 1 per field of view 
on average. This compares to a probability of 0.3 for 
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finding a single one of the se objects in the cu rrent largest 
area HST survey (BoRG,|Trenti et aL|(|2011|)). The rota¬ 
tion curves will be exciting targets tor I'hirty Meter class 
ground based telescopes. Forming from lightly enriched 
material, these primeval disk galaxies are special objects. 
Observing them and comparing them with simulations 
and analytic theory should enable us to make inferences 
about dark matter, the role of angular momentum and 
the fundamental principles of galaxy formation (such as 
tidal torques, and smooth accretion of unenriched mate¬ 
rial) which are not possible in the later, more complex 
Universe. 


We acknowledge funding from NSF OCI-0749212 and 
NSF AST-1009781. The BlueTides simulation was run 
on Blue Waters supercomputer at the National Center for 
Supercomputing Applications. 
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